Continuum Hartree-Fock-Bogoliubov theory for weakly bound deformed nuclei using 

coordinate-space Green's function method 
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We formulate a new scheme of the Hartree-Fock-Bogoliubov mean-field theory applicable to weakly 
bound and pair correlated deformed nuclei using the coordinate-space Green's function technique. 
On the basis of a coupled-channel representation of the quasiparticle wave function expanded in 
terms of the partial waves, we impose the correct boundary condition of the asymptotically out- 
going waves on the continuum quasiparticle states. We perform numerical analysis for '^^Mg to 
illustrate properties of the continuum quasiparticle states and the pair correlation in deformed 
nuclei near the neutron drip-line. 
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I. INTRODUCTION 



The Rl-beam facilities in the new generation will en- 
large significantly the experimentally accessible region in 
the nuclear chart, in particular in medium and heavy 
mass domains. An interesting area may be the 10 ^ Z < 
20 and N>2Q region, where the selfconsistent mean- 
field theories predict that the shape deformation sys- 
tematically occurs even if nuclei are close to the drip- 
line P, 0, y, 0, H, @, 0| ■ It is the presence of weakly bound 
neutrons that makes these nuclei interesting, and further- 
more the possible shape deformation will bring about ad- 
ditional mechanism influencing the single-particle motion 
and the many-body correlations such as the pairing and 
the collective excitations. 

A promising theoretical framework to describe this sit- 
uation may be the selfconsistent mean-field approachjH, 
1^, ■ We consider more specifically the Hartree-Fock- 
Bogoliubov (HFB) method to describe the pair correlated 
and deformed ground state [Sj, llO|, and the quasiparticle 
random phase approximation (QRPA) to describe exci- 
tation modes built on the ground state[^, Note here 
that, since the nucleons are bound only weakly and the 
threshold energy for the nucleon separation is low, one 
has to formulate the HFB and QRPA methods on the 
basis of correct description of the asymptotic forms of 
the wave functions of weakly bound and unbound con- 
tinuum quasiparticle states [ll|, [H, [H, [13, [3 ■ The for- 
malisms that fulfill this requirement, which we shall call 
the continuum HFB[i1 [ll, El [H or the contin- 
uum QRPAjl4l. lisj . are limited mostly to spherical nu- 
clei. Therefore we need to extend, as the first step, the 
continuum HFB to deformed nuclei, i.e., we need to for- 
mulate the deformed continuum HFB theory. Hamamoto 
(20I [2l| has analyzed in detail the quasiparticle motion 
in deformed Woods-Saxon potential by solving the HFB 
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equation in the coupled-channel representation imposing 
the boundary condition of the correct asymptotics. Re- 
cently Stoitsov et al. have introduced a formulation 
of the deformed continuum HFB utilizing the Poschel- 
Teller-Ginochio basis. In the present paper, we intend 
to give another new formulation of the deformed con- 
tinuum HFB method by extending the Green's function 
approach [isj. Our eventual target is not just to formu- 
late the deformed continuum HFB, but also to formu- 
late a continuum QRPA for deformed nuclei. Although 
the latter is not pursued here, the deformed continuum 
QRPA can be easily formulated once the quasiparticle 
Green's function for deformed nuclei is constructed in 
the coordinate representation J^J. 

The kernel of the present formulation is the quasipar- 
ticle Green's function (called also the HFB Green's func- 
tion) which satisfies the correct boundary conditions of 
the continuum quasiparticle states. This can be achieved 
by utilizing the coupled-channel representation [2^ |2l[ of 
the quasiparticle Schrodinger equation (the HFB equa- 
tion) based on the partial wave expansion. The exact 
form of the HFB Green's function is known for spher- 
ical nuclei 13] where the channels - the partial waves - 
decouple, but what we need is the one in deformed nu- 
clei for which there exists the coupling among the par- 
tial waves. We mention here that the Green's func- 
tion for deformed potentials has been utilized in describ- 
ing other physical systems, e.g. the electronic response 
in molecules [2^ [23| and electrons in matter scattering 
on deformed ion potentials [13], for which the coupled- 
channel representation based on the partial wave expan- 
sion is also employed. Indeed the exact form of the 
Green's function in the general coupled-channel system 
[13] can be extended to our problem, i.e. to describe the 
quasiparticle wave functions in pair correlated deformed 
nuclei. In the HFB theory, one also needs to calculate the 
density and the pair density, or the generalized density 
matrix in general, by summing up the wave functions of 
all the quasiparticle states including those in the contin- 
uum. The Green's function formalism fis'] can be utilized 
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also to efficiently perform this summation. In this way we 
obtain a complete scheme to calculate the HFB ground 
state and the single-particle properties influenced by the 
pair correlation. We explain the details of this formalism 
in Section II. 

As a demonstration of the deformed continuum HFB 
method, we perform in section III numerical analysis for 
a neutron-rich nucleus ■^^Mg, which is situated near the 
neutron drip-line and is predicted to be prolately de- 
formed in many selfconsistent mean- field calculations [l|, 
H, H, 0, H, H, I3l- The purpose of this analysis is two 
fold. The first is to investigate how the neutron pair 
correlation changes as the binding of neutrons becomes 
weaker and stronger. Secondly we would like to reveal 
peculiar properties of the single-particle motion in the 
pair correlated deformed nuclei near the neutron drip- 
line, where the coupling of the quasiparticle states to 
the continuum orbits may introduce new features in the 
single-particle motion. We shall compare our results with 
those of Hamamoto [20|, IIH , who performed a pioneering 
analysis of the continuum quasiparticle states in a for- 
malism satisfying the correct asymptotic forms, but not 
on the basis of the selfconsistent treatment of the pair 
correlation. Section IV is devoted to the conclusions. 



II. DEFORMED CONTINUUM HFB THEORY 
IN THE GREEN'S FUNCTION FORMALISM 

A. Coordinate-space HFB equation 



It obeys the HFB equation 



h-X 
h 



h 

-h + X 



0(1) (ro-, E) 
(t)^^^rcr, E) 



E 



(j)^^\ra,E) 
0(2) (rCT,^;) 



(2) 

where E is the quasiparticle excitation energy, and A is 
the chemical potential or the Fermi energy which should 
be determined to constrain the expectation value of the 
nucleon number. We omit the isospin index for simplic- 
ity of notation. Here and hereafter, we assume the time 
reversal invariance of the HFB ground state \^), the as- 
sociated Hartree-Fock Hamiltonian h and the pair Hamil- 
tonian h. We also make a simplification that the Hartree- 
Fock potential and the pair Hamiltonian are local ones: 



h = ^ — V^ + V{rs) 
2m 



h = Air), 



(3) 



as is realized in the cases of some Skyrme effective inter- 
actions, e.g. SkP[l^ with the effective mass m* = 1, or of 
the phenomenological Woods-Saxon potential. The pair 
Hamiltonian becomes a local pair potential A(r) when 
the contact force is adopted for the pairing interaction. 
The position dependent effective mass m*{r) is often en- 
countered in the Skyrme Hartree-Fock(-Bogoliubov) the- 
ories. An extension to this case is straightforward, but 
we do not deal with here. 

We expand the quasiparticle wave function with re- 
spect to the partial waves specified by the angular quan- 
tum numbers jlm, abbreviated by L hereafter: 



The Bogoliubov's quasiparticle plays a central role in 
the HFB theory [13]. The wave function of the quasi- 
particle state has two components, and is written in the 
coordinate-space representation as 



^^\ra,E) 
'•^\ra,E) 



(t>^l\r,E)YL{fa) 
4>f{r,E)YL{ra) 



(4) 



where YL^fa) is the spin spherical harmonics. The 
HFB equation is then transformed to a coupled-channel 



(1) formj2d. l2lj for the radial wave functions (p)^ {r,E) 



2 d i{e + 1) 



2m dr^ 2m r dr 2m 



-X-E 



2m dr^ 2m r dr 2m 



2 d _^h' £{£ + 1] 



L' 



(5) 



+ J2ALL'ir)AV{r,E) = 0. 



(6) 



Here the "channels" are labeled by the quantum number ufj^,{r) , u^lj^,{r) and ALL'{r), which are defined by 
jlm{— L) and the index i = 1,2 specifying the upper 

and lower components of the quasiparticle wave func- "li,/ ('') + w^"'^, (r) = / df'S^Yl[fa')V{rs)Yi^:{fa) 

tion. The coupling among the channels is governed by J „„, 



(7) 
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^LL'{r) = I dfY,Yl{fa)A{r)YL'{fa). 



(8) 



Note that the coupled-channel representation is often em- 
ployed to describe scattering states in the non-spherical 
potential problems [H [H [H ill . We have to trun- 
cate the partial wave expansion in practical calculations, 
and we denote N to represent the number of partial waves 



J 



to be included. 

Using the two component radial wave function 



the coupled-channel equation can be written as 



(9) 



^Mr,E) + ^^cbLir,E)+Y^v'il,ir)^cbL'{r,E)+Y^v^l!l = 0, 

where we introduced 2x2 matrices 



v^Lkr,E) = 



au^^j^, (r) + aX + aE 
aALL'{r) 



-aALL'ir) 
e(e+i) (0) / \ I \ 77 



(10) 



(11) 



(1) / \ / -au\i,{r) 

\ -au^i,,(r) 



(12) 



and a constant a = ImjTi' . It is also possible to write 
this equation in a grand matrix form 

|l+2|.+^(i)(,)|. + ^(0) E)\<i>{r,E)^^, 
dr r or or J 

(13) 

where the radial wave functions with different L's are 
combined to form a 2N dimensional vector 



{r,E) 



( Mr,E) \ 
<l>LAr,E) 



(14) 



\<j>LAr,E) J 
and 2N x 2N matrices v^"^ and v^^^ are defined by 

(n=l,2) (15) 



B. Boundary conditions for the quasiparticle wave 
function 

The quasiparticle states with the energy higher 
than the nucleon separation energy form a continuum 
spectrum[ll|,|ll,[ll|. They are the states with \E\ > |A| 
while the discrete quasiparticle states lie in the energy 
range \E\ < \X\. In order to describe the continuum 
quasiparticle states, it is an essential condition that the 
quasiparticle wave functions must have correct asymp- 
totic forms at far outside r oo where the potentials 
vanish. 

We therefore impose the boundary condition at r — > oo 
that the radial wave functions should be connected to the 



following asymptotic form. There are 2N independent 
solutions where the upper or lower component of a given 
partial wave L' is dominant: 



Mr,E) 



r 





'LL' 







Here 
iJ+(fc,r) = 



(16) 



ji(fc, r) -f ini{k, r) (for neutron ) 
Fi^z{k,r) + iGi^z{k,r) (for proton ) 

(17) 

is the Hankel/Coulomb function. The wave numbers 
in Eq. ^ are given by k±{E) = yj2m{\±E)/h and 
their branch cuts are chosen so that lmk±{E) > is 
satisfied; 13]. This boundary condition is equivalent to 
imposing that the wave functions are connected to the 
asymptotic out-going or exponentially decaying waves. 
Combining these solutions in columns, we introduce a 
2x2 form tp^^^ (r, E) satisfying 



(out), 



i^±fl5,„ 

^^^5,,, 



(18) 

and similarly a 2N x 2N matrix form <3>(°"*)(r, E) defined 

by 



[^^'^-'\r,E)U'^^Tr'ir,E) 



(out). 



(19) 



At the center of the nucleus, we consider the radial 
wave functions which arc regular at r = 0: 



(in) I 



^0 ( r'5LU 



(20) 



roLu 

We denote it also in the 2N x 27V matrix form as 

[^'^"^\r,E)]LL'^^TL'ir.E). (21) 
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C. HFB Green's function in the coupled-channel 
representation 

Now we introduce the Green's function defined for the 
coordinate-space HFB equation It is expressed for- 
mally by 



E- 



h-X h 
h -h + X 



(22) 



and may be denoted 

(23) 

in the coordinate-space representation. Note that the 
HFB Green's function G{E) has the 2x2 matrix form, 
whose diagonal and off-diagonal components are often 
referred to as the normal and abnormal Green's function, 
respectively. Using the partial wave expansion it may be 
expanded as 



G{ra,r'a',E) ^J2^L(^'^)gLL'{r,r',E)Y£,{f'a'), 



LL' 



We can also introduce the 2N x 2N matrix form g(r, r', E) 
of QLL'ir, r',E): 



(24) 



.(12)/ 



[g{r,r',E)]LL' =9LL'{r,r',E). 



(26) 



It is easy to derive from the definition of the Green's 
function that g{r,r',E) must satisfy the equation 



1' 
dr 



U -I- + v(^) (r)|- + v(") (r, E)] g(r, r', E) 
r or or J 

6{r — r'' 



where J is defined by 
J = 



1 

-1 



aJ, (27) 



(28) 



with Im being the N dimensional unit matrix. 

We now seek the Green's function g(r, r', E) which sat- 
isfies both Ea. (p7)l and the boundary conditions, ([T8|) and 
pop , imposed at r — > 00 and r = 0, respectively. Accord- 
ing to Rcf. 24] , we write it in a form 



gLL'{ry,E) = (ctrl{r',E)^fi^„{r,E)e{r' - r) 

L" 

+ (/, E)^^ll^} (r, E)e{r~ r')) (,29) 



or equivalently 



;(r, r',E)^ (r, £;)C(™'^(r', E)ey - r) 

-h$(°"*)(r, £;)C(°"*)^(/, E)e{r ~ /), (30) 



where C(™^(r, i?) and C(°"*)(r, i?) are matrix functions 



1 1 out) 



yE)] 



^{in/out) 



LL' 



yE). 



(31) 



Substitution of Eg. ([50]) into Eq. (P7)l leads to the equa- 
tions for C(™)(^£^) and y°''*^r,E): 



(Jiout) ^^^^^ 



3 m 



M[E) 



which corresponds to Eq.(23) of Ref.|2j|. If the Hamil- 
tonian does not contain the coupling term with the first 
derivative, i.e. v(i^(r) = 0, the solution has a simple form 

C(")^(r,£;) = a(M^)"^$(°"*)^(r,£;), 
Cio^y^r.E) = aM-i$(")^(r,£;), (33) 
with the Wronskian matrix defined by 

$(™)'^(r, E)J (^-^$(°"*)(r, E)^ 

A$(™)^(r^£;) ) j$(°"*)(j.,£:) , (34) 
dr J \ 

and consequently the Green's function is written as 

g{r,r',E) ^ a (r, S) (M^)"^$(°"*)^(r',£;)6l(r' -r) 

+ $(°"*)(r,£;)M-i$(™)^(r',£:)6'(r-r')) , (35) 

in parallel to Eq. (36) of Ref.^. 

In the following numerical application, however, we do 
not use Eq. fSS)) . but instead we solve directly Eq. ([5^ to 
obtain C^^H?": E) and C(°"*)(^, E). We then use Eq. (i29l) 
to evaluate the HFB Green's function. This procedure is 
applicable to the case where the coupling with the first 
derivative is present. We found also that the procedure 
provides us better numerical accuracy than to use Eq. |35)) 
even when there is no first derivative term. 

D. Generalized density matrix 

Key quantities in the HFB theory are the normal den- 
sity matrix 

p{ra,r'a') = (1'|^^(r'cr')V^(rcr)|*), (36) 
and the abnormal density (pair density) matrix 

p{ra,r'a') = {^\i}y d')i){ra)\^) , (37) 
which can be combined as the generalized density matrix 

^ ' ' \p*{rc7,r'a') Srv'S^^' - p{r(T,r'a') J 

(38) 



Let i be the index to specify the quasiparticle states, 
and (j)i{ra) be the eigen solutions of the HFB equation 
(l2|). Then the generahzed density matrix is obtained by 
summing up products of the quasiparticle wave functions 
over all the quasiparticle states ,141: 

R{ra, r'a') = ^ 4>lira)^{r'a') (39) 

i 

where is a conjugate wave function of 0^, defined by 

Here a implies (f>^l\r(7) = {—2a)if)[^\r~a), and (f)~^\ra) 
is the time-reversal (j)~.^\ra) = T(j)f\rfj) — (jyf'* {ra) of 
</)«(ra). 

We perform this summation by using a contour integral 
of the HFB Green's function in the complex plane of the 
quasiparticle energy E flS*! . Then the generalized density 
matrix is given by 

i?(rcr,rV) = — / G(rcr,r'a\E)dE. (40) 

We use the partial wave expansion also for the generalized 
density matrix: 

i?(ra,rV) = J2YL{fa)RLL'{r,r')Yl,{f'a'), (41) 

LL' 

RLL'{r,r')^^ I gLL'ir,r',E)dE. (42) 

The contour C should enclose the negative energy side of 
the real E axis. The use of the contour integral enables us 
to impose the proper boundary condition for the contin- 
uum quasiparticle states using the HFB Green's function 
gLL' described in the previous section. 

In order to make the numerical integral efficiently, we 
choose C as a circular path shown in Fig[T] although a 
rectangular contour is adopted in Ref.[l3|. Then the in- 
tegral reads 

RLL'{r,r') = ^ l'\LL'{r,r',E{OW^^-^^dC, 

(43) 

where Ecut is the maximal quasiparticle energy which 
define the cut-off of the sum \E\ < Ecut of the quasipar- 
ticle states. The circular path makes the integrand to 
behave smoothly. In actual numerical calculations, we 
split C into two semicircles < ^ < tt and tt < C < 2tt, 
and then apply the higher-order Gauss-Legendre quadra- 
ture with points to perform the numerical integra- 
tion in each semicircle {Mq points in total). The local 
density p{r) = Pi''''^^ ''''^) the local pair density 
p{r) = J2a 'f''^) calculated accordingly. 
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lm{E) 




FIG. 1: The circular path C adopted to perform the contour 
integral (|43p . The thick solid line represents the negative 
energy continuum quasiparticle states with E < X while the 
crosses represent the negative energy discrete quasiparticle 
states X < Ei < 0. 



Given the methods to calculate the quasiparticle wave 
functions and the generalized density matrices, an usual 
iteration procedure can be applied to obtain a converged 
HFB ground state I^P), the associated quasiparticle states 
and the densities. The Fermi energy A is also determined 
selfconsistently so that the expectation value of the neu- 
tron/proton number is constrained to N/Z of the nucleus 
under consideration. 



III. APPLICATION TO A DEFORMED 
DRIP-LINE NUCLEUS 

In the following we shall illustrate the deformed contin- 
uum HFB method with numerical examples. We would 
like to discuss some features of the quasiparticle spec- 
tra and the pair correlation which are characteristic to 
deformed nuclei near the neutron drip-line. As an ex- 
ample, we choose ^^Mg which is predicted to be pro- 
lately deformed and close to the neutron drip-line in the 
mean- field calculations [Tl, H Uli, H, H, 0, 113 ■ Experimen- 
tally, ■^^Mg is unbound and "^°Mg is the most neutron-rich 
bound isotope identified so far [2^. We do not deal with 
^"Mg for which the prediction of the deformation is more 
subtle [l, i2|, i3|, 1^, ,5„ ^, i?,, i27,] . 

A. Model and numerical procedure 

In the following analysis we intend to clarify qualitative 
features rather than to make precise quantitative predic- 
tions of the specific nucleus. We simplify the Hartree- 
Fock potential by replacing it with an axially-symmetric 
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deformed Woods-Saxon potential 



(44) 



where f{r,e) = (i + e-('-«(9))/io)-i and R{e) ^ Ro{l + 
/3Y2o{0)), together with the spin-orbit potential 



r dr 



£-s. 



(45) 



We choose the deformation parameter P = 0.3, being a 
typical value of the mean- field predictions [U, [|, i, i, d, 
a, 0, [13] ■ The other Woods-Saxon parameters are taken 
from Bohr and Mottelson[2^. We neglect the deforma- 
tion in the spin-orbit potential. The neutron Nilsson di- 
agram of the deformed Woods-Saxon potential is shown 
in FiglJl If there is no pair correlation, the highest occu- 
pied orbit is the one labeled with the Nilsson asymptotic 
quantum number [310] ^ 



<D 

o 
t: 

ni 
a. 

cn 




-2 
-4 
-6 
-8 
-10 
-12 
-14 



[312] f 



1 


-—-^ I Z Z = = t = - - _ J 




[20 


] 




~~ -~ — 2_ 






Si 




2 






38Mg 8 

1 1 


1 1 





1310' 



-0.2 -0.1 



0.1 



0.2 



0.3 0.4 



FIG. 2: The neutron Nilsson diagram of the adopted Woods- 
Saxon potential for ^*Mg. The neutron Fermi energy is lo- 
cated near the orbits labeled with the Nilsson asymptotic 
quantum numbers [310] | and [312] |. The other orbits with 
f2 = 1/2 are also marked with the arrows. 



Concerning the pairing properties, we derive it in a 
selfconsistent way from an effective pairing interaction, 
for which we adopt the density-dependent delta interac- 
tion (DDDI), 



The selfconsistent pair potential is then given by 
A(r)4.o(l-.(^)")p(.). 



(46) 



(47) 



Here we adopt the parameters {vq = —458.4 MeVfm ^, 
r, = 0.76, a = 0.59, pc = 0.08 fm"^) given in Ref . [M, HH , 
where vq is determined so that the DDDI describes the 



scattering length a = —18.5 fm under the given cut-off 
energy Ecut- 

The selfconsistent pair potential depends on the pair 
density p{r), and it varies in different physical situations. 
To provide a contrast to this case, we perform another 
calculation where the pair potential is replaced by a fixed 
phenomenological one having a Woods-Saxon shape 



(48) 



as adopted in Refs. |20l. l2lj . According to Refs.|20l. 
[2l| we constrain the depth of the pair potential 
by the value of an average gap A defined by A = 
JdrA{r)fir,e)/Jdrfir,e). 




10 12 14 16 18 20 22 24 26 
2nmax 

FIG. 3: Dependence of the neutron pair interaction energy 
Epair on the maximal value Qmax of the magnetic quantum 
number of the quasiparticle states. 



The numerical details are as follows. We consider only 
neutrons. Since we assume the axial symmetry, the pro- 
jection of the total angular momentum along the z-axis 
is conserved. The quasiparticle states have the magnetic 
quantum number Q , the eigenvalue of jz, and we solve 
the HFB equation separately for each il. The coupling 
potentials u^^j^,{r) and A^j^j^,{r) are evaluated according 
to Eqs.([7|) and ([5]), where the integral over the azimuthal 
angle 9 is performed numerically by means of the higher- 
order Gauss-Legendre quadruture with Ng points. The 
coupled-channel HFB equation (fT^]) is solved by means 
of the Runge-Kutta-Nystom method [33| in the region 
r = — Rmax with Rmax = 15 fm and an equidistant 
^ interval Ar = 0.2 fm. We impose the boundary condi- 
' lion that the wave function is connected to the asymp- 
totic form, Eg. p^ . at r = Rmax- We use the code 
cwf complex 33] to calculate the Hankel and Coulomb 
functions at complex energies. The cut-off quasiparticle 
energy is Ecut — 60 MeV. We also truncate the quasi- 
particle states in terms of the magnetic number by in- 
cluding those with —Qmax < < ^max- Figure [3] shows 
the dependence of the neutron pair interaction energy 
Epair = 5 / drA{r)p{r) on the choice of flmax- A rea- 
sonable convergence is obtained with rimaa;^ 21/2, and 
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hence we adopt il.max = 21/2 in the following analysis. 
The parameters for the Gauss-Legendre quadrature are 
Ne = 60 and A/<^ = 100. 



B. Density profiles 

The deformed continuum HFB method can describe 
the exponentially decaying asymptotics of the densities. 
Figure [4] shows the density profiles of the HFB ground 
state. Here we make the multipole expansion of the den- 
sities as 



p{r) 



p{r) 



- ^p,(r)y,o(r), 
= ^p,(r)y,o(r), 



(49) 



(50) 



and plot the monopole and quadrupole parts 1 — 0,2. 
For comparison, we also show the results of another cal- 
culation where we impose the box boundary condition 
4>L{Rmax) = at the boundary r = Rmax- In the lat- 
ter case the whole quasiparticle spectrum including the 
continuum part is discretized. It is clear that the calcu- 
lation with the correct boundary condition can describe 
the exponentially decaying asymptotics both for the den- 
sity and the pair densities while the calculation with the 
box boundary condition fails. 

We remark here that both the monopole and 
quadrupole parts exhibit the same slope in the expo- 
nentially decaying tail, i.e., the ratios P2{'r)/po{f) and 
P2 (r) / po (r) are approximately independent of r in the 
tail region. This means that we can define the "deforma- 
tion" of the exponentially decaying tails in terms of the 
shape of the equi-density contour curves. Note however 
that the size of the deformation in the tail differs from 
that of the nuclear surface ( which may be defined by 
the deformation of the equi-density contour at the half 
central density). It is also seen that the tail deformation 
in the density p{r) and that in the pair density p{r) are 
different. It is interesting to quantify the deformations in 
the tail region, but we leave it for future works since in 
the present calculation we do not make the selfconsistent 
treatment of the deformation of the Hartree-Fock poten- 
tial. The calculation in Ref.[l3] indicates also that the 
exponentially decaying tail has the same slope in different 
directions. 



C. Weak-binding eff'ect on pair correlation 

We can regard the calculated -^^Mg as a quite weakly 
bound system since the neutron Fermi energy in the 
present calculation is just A = —0.889 MeV. It is in- 
teresting to see what happens if the system is weakly 
bound further more, or if the system is more strongly 
bound. To see the effects of the weak binding, we inves- 
tigate the pairing correlation by varying artificially the 
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FIG. 4: (Top) The monopole and quadrupole parts po{r) and 
P2{r) of the neutron density p{r). The results of the calcula- 
tion using the box boundary condition are also plotted with 
the dashed curve for comparison. (Bottom) The same but for 
the monopole and quadrupole parts po{r) and P2{r) of the 
neutron pair density p{r). 



depth parameter V^^ of the Woods-Saxon potential as 

+ « (51) 

where we choose a = +1, 0, — 1, — 2, — 3 and 
—4 MeV. The neutron Fermi energy is A = 
-0.471,-0.889,-1.335,-1.802,-2.289 and -2.793 
MeV for a ^ +1,0,-1,-2,-3 and -4 MeV, respec- 
tively. 

Figure [H] shows the monopole parts po(r) and po{r) of 
the density and pair density, respectively, of neutrons for 
different values of a. It is seen that the exponentially 
decaying tails of the density and the pair density ex- 
tends (shrinks) as the Fermi energy becomes shallower 
(deeper). This is naturally expected as the exponen- 
tial decay behaves kinematically as p{r) oc e"^'"' and 
p(r) cx e~'^^ with the exp onent k — y/2m\X\/h related 
to the Fermi energy X\M- We also see another notice- 
able feature that as A becomes shallower the magnitude 
of the pair density increases not only in the tail region 
but also in the whole region of r, including the interior 
and the surface areas. Apparently the pair correlation 
is enhanced by the weak binding. To characterize the 
overall magnitude of the pair correlation we evaluate the 
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FIG. 5: (Top) The monopole part of the neutron density po{r) 
for various shifts of the Woods-Saxon potential depth with 
V°s + a [a = +1,0, -1, -2, -3 and -4 MeV). The inset is 
the same plot but in the log scale. (Bottom) The same but 
for the monopole part of the neutron pair density po{r). 



average pairing gaps A^y and A^v defined by 

^ ^ JdrA{r)p{r) 
I drp{r) 

^ ^ JdrA{r)p{r) 
J drp{r) ' 

and the pair interaction energy 



drA{r)p{r) 



(52) 



(53) 



(54) 



They are plotted in Fig[6l where the enhancement due to 
the weak binding is seen directly in these quantities. 

In order to investigate the origin of the enhanced pair- 
ing in the weak binding cases, we compare with another 
calculation where the pair potential is replaced by a phe- 
nomenological one, Eq. (P5|) . with fixed strength A = 1.5 
MeV. The calculated average pairing gaps and Auv 
and the pair interaction energy Epair are shown also in 
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FIG. 6: (Top) The average pairing gaps A„„ and Auv plotted 
with the solid curves as a function of the neutron Fermi en- 
ergy A for the selfconsistently calculated pair potential. The 
results for the fixed phenomenological pair potential is also 
plotted with the dashed curves. (Bottom) The neutron pair 
interaction energy Epair as a function of the neutron Fermi en- 
ergy A for the selfconsistently calculated pair potential (solid 
line) and also for the fixed phenomenological pair potential 
(dashed line). 



Figl6l where, in contrast to the selfconsistent case, A„^ 
and A„„ stay almost constant or even decreases slightly 
with lifting the Fermi energy. The monopole density 
Poif) and the monopole pair density po(r) are shown in 
FigEl It is seen that the pair density inside the nucleus 
decreases as the Fermi energy becomes shallower. This is 
opposite to the trend seen in the case of the selfconsistent 
calculation. Consequently there is very little influence of 
the weak binding on the pair correlation provided that 
we fix the pair potential. In other words, the enhance- 
ment of the pair correlation caused by the weak binding is 
described only when we calculate the selfconsistent pair 
potential. In the selfconsistent case, the increase of the 
pair density in the tail region induces the increase of the 
pair correlation not only in the tail region but also in the 
entire region. 
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FIG. 7: The same as Fig[5]but for the fixed phenomenological 
pair potential. 



D. Structure of continuum quasiparticle states 

In the HFB theory, both the neutron density p{r) and 
the expectation value (N) = J drp{r) of the neutron 
number can be expressed as summations of contributions 
of the individual quasiparticle states. Since the quasipar- 
ticle states are specified by the energy eigenvalue E, we 
can define the contributions per unit energy (the density 
per energy) to the neutron number {N) . We denote it by 
n{E), and call the occupation number density. It satisfies 

(N) = /g dEn{E) by definit ion, and can be calculated 
in terms of the HFB Green's function by 



i{E) = -ImY^ J drG^^^\r(7,r(j,-E -ie). 



(55) 



This is the level density of the quasiparticle states, but it 
is weighted with the occupation number (corresponding 
to the factor in the case of the BCS approximation). 
Similarly we can define the pair number density n[E) by 




i{E) = -Im^ j drG'^^^\rcr,ra,-E -ie) 



(56) 
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FIG. 8: The occupation number density n{E) (top panel) 
and the pair number density n{E) (bottom) for the neu- 
tron quasiparticle states with Q = 1/2. The solid curve is 
the results of the deformed continuum HFB calculation with 
the correct boundary condition on the asymptotics while the 
dashed curve is obtained by imposing the box boundary con- 
dition (f>{Rmax) = 0. The arrow indicates the threshold energy 
|A| = 889 keV for the continuum quasiparticle states. 



for which the level density is weighted with the pairing 
factor (the uw-factor in the BCS). n{E) represents the 
contribution of the quasiparticle state at energy E to the 

pair number {^-^^ — J drp{r). Note that the constant 
e plays a role of the smoothing parameter. If the quasi- 
particle spectrum contains discrete states, there emerge 
delta functions in n{E) and n{E) at the energies of the 
discrete quasiparticle states. Having the imaginary part 
ie in Eqs. (|55)) and the delta function becomes a 

Lorentzian function with FWHM of 7 = 2e. In the fol- 
lowing, we take e = 25 keV (7 — 50 keV), and we shall 
investigate the structure of the quasiparticle spectrum in 
terms of n{E) and h{E). 

Figure [5] shows n{E) and n{E) for neutron states with 
il = 1/2. Note that the threshold energy of the contin- 
uum spectrum is Eth = |A| = 889 keV. It is seen that the 
continuum spectrum above Eth consists of several res- 
onances forming narrow peaks, and a broad continuum 
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background for which we hardly assign resonance struc- 
tures. We can identify a resonance a,t E — 1.24 MeV as 
a quasiparticle state corresponding to the orbit [310]^ of 
the deformed Woods-Saxon potential (cf. Fig[21). Nar- 
row resonances located at E ~ 3.36, 4.61, 7.36 and 11.68 
MeV correspond to the other = 1/2 orbits marked with 
the arrows in FigO The widths of these resonances are 
smaller than the smoothing width 7 = 50keV except the 
lowest energy resonance at E — 1.24 MeV, which has a 
width apparently larger than those of the others. Note 
that the [310]^ orbit is a bound state if we neglect the 
pair potential, and the finite and relatively large width 
is brought by the pair potential. 

Another important feature seen in FiglB] is that the 
non-resonant continuum states have significant contribu- 
tion to the pair number density n{E) in an extent com- 
parable to that of the narrow resonance states. This 
should be contrasted to that the contributions of the 
non-resonant continuum states to the occupation num- 
ber density n{E) is very small. The small contribution 
of the non-resonant continuum states to the occupation 
number density n{E) is pointed out in Ref. 21], but we 
emphasize here that this is not the case for the pair num- 
ber density n{E). 

In Fig. [H we also show the results obtained with the 
box boundary condition 4>{Rmax) = 0. In this case the 
quasiparticle states are all discretized. It is obvious that 
the discretized quasiparticle states obtained with the box 
size Rmax — 15 fm fail to describe the width of the res- 
onances, for instance, that of the [310] |. Obviously the 
discretization is too crude to describe the continuous be- 
havior of the non-resonant continuum states. 

In Figs. [9] and [10] we show and h{E) combin- 

ing several values of Vl. We mention here two features, 
(i) Firstly, the lowest energy resonance (with Vl = 1/2) 
and the second lowest one (with = 5/2) have signif- 
icantly different widths even though the centroid ener- 
gies of both resonances are comparable. We note here 
that the single-particle energies of the orbits [310] ^ and 
[312] I , which correspond to the lowest and the second 
lowest resonances, respectively, are also very close (cf. 
Figl2]). The difference in the width may be attributed 
to the difference in the r nag netic quantum number 17 in 
a way discussed in Ref. [26|: the centrifugal barrier for 
the [310] I state is lower than that of the [312] | state 
since the former contains the partial wave with I ~ 1 
while the latter does not. (ii) Secondly, it is seen in the 
pair number density n{E) that the contribution of the 
non-resonant continuum states is larger for states with 
smaller f2. Note, however, that the relative importance 
of high/low- non-resonant states depends on the quasi- 
particle energy. The low-fi states such as Q, — 1/2 and 
3/2 are dominant as far as the non-resonant quasiparti- 
cle states at low energies ^ 10 MeV are concerned (cf. 
FigE]). Concerning the non- resonant continuum states at 
higher energies £' <^ 30 MeV, however, we see in FigfTUl 
that contributions from various il's do not depend very 
strongly on fi: there is only a small difference by a fac- 
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FIG. 9: (Top) The occupation number density n(£') of 
neutrons shown separately for different values of SI = 
1/2, 3/2, ■ ■ ■ , 11/2. The arrow indicates the threshold energy 
|A| for the continuum quasiparticle states. (Bottom) The 
same but for the pair number density n{E). 



tor of 2 in n{E) between J7 = 1/2 and 11/2 states at 
E^30 MeV. In other words, the high-fi states also have 
sizable contributions to the pair correlation. 

The decomposition of the quasiparticle states into the 
partial waves is another useful tool to analyze the struc- 
ture of the quasiparticle states. In fact we can decompose 
the occupation number density n{E) and the pair num- 
ber density n{E) into different partial waves: 

niE) = ME). n{E) = (^7) 

L L 

where the contribution from each partial wave L is given 

by 

nL{E) = -Im / drr^g'il\r,r, -E - ie) (58) 



and similarly for hL{E). The partial wave decomposition 
of the fi = 1/2 states is shown in FigfTTl It is seen that 
each of the resonances has specific partial wave contents 
reflecting the character of the resonances. For instance, 
the lowest energy resonance [310] 5 consists of dominant 
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FIG. 10: The same as Fig|9]but in the log scale and in a wider 
energy range. 



Pi/2 ~ P3/2 components and sub-dominant /7/2 compo- 
nent. The second resonance consists of dominant fj/2 
and sub-dominant ^3/2- 

We can also argue the structure of the non-resonant 
quasiparticle states in terms of the partial wave contents 
in h{E). For instance, the continuum states in the inter- 
val 1.5 ^i?^ 3.0 MeV is characterized by the relatively 
largest component of P3/2 waves and sub-dominant 
/7/2 ~ /5/2 components. We see also some structures in 
the intervals 2.5 ^£'^5.5 MeV and 5.5 <£'^ 8, which 
are characterized the relatively largest components /5/2 
and 59/2, respectively. These structures presumably orig- 
inate from resonance states in the deformed Woods- 
Saxon potential lying high above the Fermi energy. They 
coexist in the same energy region with smooth (almost 
flat in E) components containing several partial waves 
with I — Q ~ 4. At higher energies, say E^30 MeV, 
there remain only smooth non-resonant components in- 
cluding all the partial waves up to very high orbital an- 
gular momentum 1^6 — 10. The contributions of high-Z 
orbits is consistent with the high-fJ components observed 
in FigHni 
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FIG. 11: The partial wave decomposition of the occupation 
number density n{E) (top panel) and the pair number density 
n{E) (bottom) of neutron Q — 1/2 states. Each curve labeled 
with Ij represents the partial sum X^^^^j^. riLiE) (fiL^E)) con- 
taining from Si/2,Pi/2,P3/2, • ■ • to Ij. 



E. Weak-binding effects on quasiparticle motion 

In this subsection we shall elucidate how the weak 
binding influences the neutron quasiparticle motion and 
its coupling to the pair correlation. For this purpose 
we first analyze how the neutron pair number density 
n{E) evolves when we vary the Fermi energy by shift- 
ing the depth V^g of the Woods-Saxon potential as 

^ ws ^ ws ' ^ 

Figure [T^] shows the pair number density fi{E) of 
— 1/2 and 5/2 states calculated for various Woods- 
Saxon potential depths. It is seen first of all that the 
influence of the weak binding is stronger on the 0=1/2 
states than on the = 5/2 states. (This is along with the 
observation in Refs. [13, [lH, [2^.) The influence is espe- 
cially strong on the low-lying quasiparticle states. Here 
we concentrate on O = 1/2, and focus on the lowest en- 
ergy quasiparticle state (or the quasiparticle resonance) 



12 




E + \n [MeV] 



FIG. 12: (Bottom) The pair number density h{E) of the 
neutron Q. = 1/2 states is shown for various shifts a — 
-4,-3,-2,-1,0 and +1 MeV of the Woods-Saxon poten- 
tial depth V^g — > V^g + a. The horizontal axis represents 
the quasiparticle energy \ + E whose origin corresponds to 
the threshold energy Eth = |A| for the continuum spectrum. 
(Top) The same but for the = 5/2 states. 



which corresponds to the [310]i orbit, the second lowest 
one that is visible only for a = —4, —3 and —2 MeV (and 
may be also a = —1) cases, and the non-resonant struc- 
tures around these resonance states. Note that there are 
little influences of the weak binding on the third, fourth 
and fifth lowest peaks, which correspond to the hole or- 
bits deeper than the Fermi energy by 3 ~ 8 MeV (cf. 
FigUl). 

Let us look into the [310] | state. The binding energy 
of the [310] i orbit of the deformed Woods-Saxon poten- 
tial becomes very small in the cases of a = and +1 
MeV (cf. Fig fT5|) . In these cases the quasiparticle peak 
corresponding to the [310]^ orbit is located in the con- 
tinuum region E + X > {E > \\\ — Eth)- It has a 
finite width due to the the coupling to the continuum. 
It is seen that the width increases quite steeply for the 
variation of the potential depth a = — 1— >0— >-|-l 
MeV. The steep change in the width is quite contrasting 
to that of the — 5/2 case, i.e. the [312] | state, for 
which there is only a small change in the width. Note 
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FIG. 13; The single-particle energy of the neutron Q = 1/2 
orbits in the deformed Woods-Saxon potential without the 
pair correlation for various shifts a — —4, —3, —2, —1,0 and 



+1 MeV of the Woods-Saxon potential depth 



The positive energy states are discretized with the box bound- 
ary condition at r = Rmax ~ 15fm, and shown within the 
parenthesis. The calculated neutron Fermi energy A is plot- 
ted also with the symbol +. A calculated for the fixed phe- 
nomenological pair potential is also plotted with the symbol 



also that a similar steep increase of the width can be 
seen in the second lowest peak, which corresponds to 
the Woods-Saxon orbit marked with the asterisk in Fig. 
[T3l The second lowest peak is a bound discrete state in 
the case of a = —4 MeV, and once it becomes unbound 
resonance for a > —3 MeV the width increases steeply 
along with the variation a = — 3^— 2— >— 1 MeV. 
With a = —1 MeV we see only a broad distribution or 
a wide resonance, whose width could be order of 1 MeV. 
With a = and -1-1 MeV, there remains only a smooth 
continuum which is hard to be identified as a resonance 
structure. 

The most important feature which we observe in Fig |12l 
is that the contribution of the low-lying non-resonant 
continuum quasiparticle states to n{E) increases signifi- 
cantly as the neutron Fermi energy becomes close to zero, 
especially in the cases a = and -1-1 MeV (corresponding 
to A = -0.89 and -0.47 MeV, respectively). If we look 
closely at the second lowest resonance discussed above, it 
keeps a sizable contribution to n{E) in the cases a = —2 
and —1 MeV where the width becomes very large. The 
contribution is still sizable even in the cases a = and 
+1 MeV where only a non- resonant continuum structure 
is visible as a remnant in the expected energy region. It 
is argued in Refs._20, 21] that weakly bound or resonant 
ft = 1/2 quasiparticle states contribute little to the pair 
correlation compared with that of the strongly bound 
states. Our calculation does not exhibit this behavior. 
Our analysis rather suggests that all of the weakly bound, 
resonant and non-resonant quasiparticle states contribute 
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significantly to the pair correlation in the case when the 
neutron Fermi energy is small. This is also supported by 
the observation in subsection llll CI that the weak binding 
enhances the net pair correlation. 
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FIG. 14: (Top) The pair number density h{E) of the neu- 
tron Q = 1/2 states for the shifted Woods-Saxon potential 
depths V^^ + a with a = -|-1,0 and -4 MeV. The 

arrow indicates the positions of the Fermi energy. This is the 
result of the calculation using the selfconsistent pair potential 
derived from the DDDI. (Bottom) The same but for the cal- 
culation assuming the fixed phenomenological pair potential. 
" Selfconsistent" " Fixed pair potential' gOgBa = 0, — 4gB 

Let us then investigate how the weakly bound and res- 
onant quasiparticle states couple with the pair correla- 
tion. For this purpose we shall analyze the pairing gap 
associated with these orbits. As far as the [310] | state is 
concerned, we can evaluate an effective pairing gap of this 
state in a simple way. This is because, as seen in FiglTSl 
the neutron Fermi energy A coincides within 0.1 MeV 
with the deformed Woods-Saxon single-particle energy 
Ci of the [310] I orbit for all the values of a = —4, • • • , 0, 1 
MeV. In this situation, the quasiparticle energy Ei of 
this state can be regarded^ as an effective pairing gap 



Aeff) 



since the approximate relation 



TABLE I: The peak energy Epeak of the discrete quasipar- 
ticle states and resonances in the pair number density n{E) 
for n = 1/2, and the effective pairing gaps Auv{Epeak) and 
Avv{Epeak) of the corresponding quasiparticle states, calcu- 
lated for the shifted Woods-Saxon potential depths V^^ + ct 
with a = -1-1, 0, -4 MeV. The unit is MeV for all the quanti- 
ties. In the upper half is listed the results obtained with the 
selfconsistent pair potential derived from the DDDI while the 
lower half is for those obtained with the fixed phenomenolog- 
ical pair potential, Ea. (|48|l . 
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valid in the BCS approximation leads to Ei « A^-'^''^''^'' for 

[e, - A[ < A'f^^^ We read from Fig. [Hand Table U 
that the peak position of the [310]^ resonance stays at 

Ei 1.3 MeV, and hence A^-''^^ w 1.3 MeV in all the 
cases. Namely we see that the effective pairing gap is 
almost unchanged even when — *■ 0. 

We can evaluate more directly effective paring gap as- 
sociated with the quasiparticle states. Note here that 
the density p{r, E) and the pair density p{r, E) associ- 
ated with the quasiparticle state at E is given by 



pir,E) 
p{r,E) 



ilm^ j drG^^^\ra, ra, -E ~ ie), {60} 
drG'^^'^^ra, ra, -E - ie).{61) 



= -Im 



Utilizing these quantities, we can evaluate the effective 
paring gaps of this state by 



AuviE) 



J drA{r)p{r,E) 
J drp{r, E) 

J drA{r)p{r,E) 
J drp{r, E) 



(62) 
(63) 



(59) 



The calculated values of the effective pairing gaps 
Auv {E) and A„„ {E) are listed in Table [J The val- 
ues of Auv{E) and Ayv{E) of the lowest energy reso- 
nance [310] i is consistent with the above simple esti- 
mate Af^^'' « 1.3 MeV, and we confirm again the weak 
binding does not cause significant change in the effective 
pairing gap of the [310]^ state. 



14 



In order to clarify the behavior of the effective pair- 
ing gap of the weakly bound or resonant [310] ^ state, 
we compare again with the calculations where the pair 
potential is replaced by the fixed phenomenological one 
(cf. Eg. ((48)) ). The result is shown in the bottom panel 
of Fig. [m The peak energy of the lowest energy quasi- 
particle state [310]^ varies as the Woods-Saxon potential 
depth is shifted, and we see that the effective pair- 
ing gap A^'^''^'^'' decreases as the Woods-Saxon potential 
depth becomes shallower. The reduction of the pairing 
gap of the [310]^ state caused by the weak binding is seen 
more explicitly in A„„(£') and Ayi,{E) listed in Table [H 
This reduction of the effective pairing gap due to the 
weak binding effect has been pointed out in Refs. .20, 21], 
where the mechanism of the reduction is ascribed to a 
sort of decoupling, i.e, the mechanism that the spatially 
extended wave function of the weakly bound or resonant 
quasiparticle state has less overlap with the pair poten- 
tial. But our analysis indicates that the reduction is a 
consequence of the use of the fixed phenomenological pair 
potential. The reduction of the effective pairing gap docs 
not show up if we describe the pair correlation using the 
selfconsistent pair potential. This is because the weak 
binding has two effects, one causing the enhancement of 
the pair correlation and the other causing the decoupling 
of the wave function, which have a tendency to compen- 
sate each other. 



cle wave functions of the weakly bound and continuum 
quasiparticle states. Consequently we can describe the 
quasiparticle states above the neutron separation energy 
without energy discretization. 

We have shown several numerical examples to illustrate 
effects of the weak binding and the continuum coupling 
on the pair correlation and the quasiparticle spectrum of 
neutrons in deformed nuclei near the neutron drip-line. 
The calculations are performed for ^^Mg which is cho- 
sen as an example of prolately deformed nuclei. It is 
found that there arises a significant contribution to the 
pair correlation from the non-resonant part of the contin- 
uum quasiparticle states, and the contribution grows as 
the neutron Fermi energy approaches zero. This trend is 
most significant for the quasiparticle states with f2 = 1/2. 
We confirm the stron g ef fects of the continuum coupling 
on the = 1/2 states|2ll [26j causing a steep increase of 
the width of the resonances, but we found that the con- 
tribution to the pair correlation never reduces even with 
the large continuum coupling effects. We also found that 
the effective pairing gap of the weakly bound or resonant 

= 1/2 orbits stay at sizable values. The effect of the 
decoupling from the pair correlation [l^, [2l| is not large, 
and it may be compensated by the enhancement of the 
pair correlation due to the weak binding. The use of 
the selfconsistent pair potential is crucial to describe this 
feature. 



IV. CONCLUSIONS 

We have given a new formulation of the deformed con- 
tinuum HFB method, which can be applied to deformed 
nuclei near the drip-lines. The kernel of this formula- 
tion is the use of the exact quasiparticle Green's function 
constructed on the basis of the coupled-channel repre- 
sentation in the partial wave expansion. This enables 
us to impose the correct asymptotics on the quasiparti- 
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